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I. INTRODUCTION 


ive OBJECTIVE 

The objective of this thesis is to develop a real time, 
high-resolution method for determining bottom sediment 
acoustic properties in a shallow water environment. This 
study sought to determine those properties through the 
analysis of the slope characteristics of beam reverberation 
data and higher order statistics (HOS) related to the RL 
beam time series. 

Current fleet anti-submarine warfare (ASW) operations 
are limited in shallow water in part because of their 
dependence on bottom sediment type maps that are based on 
large area averages (low resolution) derived from few 
samples [Ref. 1]. Such maps provide limited or no 
information on the bottom roughness or the sediment 
composition beneath the water-sediment interface. They 
further limit the accuracy of ASW operations because of the 
high spatial variability of the sedimentary composition, 
both laterally and vertically, that is a common feature of 
shallow water continental shelf areas. 

This thesis examines the possibility of using inversion 
techniques to deduce sedimentary properties from beam 
reverberation level (RL) data. The use of RL data is 


practical not only because a large amount of data is 


available for research, but because the RL data is gathered 
by active fleet assets rather than an oceanographic research 
vessel. In this study, the RL data was gathered during the 
Littoral Warfare Advanced Development (LWAD) Exercise 99-1, 
aboard a Ticonderoga class cruiser, using the AN/SQS-53C 


tactical sonar. 


B. SHALLOW WATER ASW PROBLEM 

Since the end of the cold war, anti-submarine warfare 
(ASW) operations have moved geographically from the deep, 
open ocean environment, to the nearshore, shallow water, 
littoral environment. Cold war ASW operations had their 
primary emphasis on locating and tracking Soviet ballistic 
missile submarines (SSBN) on patrol in the open ocean. Over 
the last 10 years, significant geo-political events have 
occurred which have shifted ASW emphasis to the littoral 
environment. With the downfall of the Soviet Union came a 
significant decrease in their SSBN patrols, while 
concurrently new threats arose in such places as the Middle 
Fast and the Korean peninsula. Smaller countries looked to 
acquire diesel-powered submarines for coastal defense. 
These vessels have become relatively affordable causing 
their numbers and countries of ownership to proliferate. 
The U.S. Navy ASW doctrine has responded to this paradigm 
change by placing increased emphasis on shallow water 


operations [Ref. 2]. 


Shallow water ASW operations are significantly more 
challenging than those conducted in the deep, open ocean. 


The littoral environment includes ocean fronts and eddies, 


stronger currents, extensive, often steep topographic 
changes, increased sea life, varied bottom sediment 
distributions, and increased shipping traffic. These 


factors combine to create a challenging ASW environment. 


Cc. REVERBERATION 

One of the most significant factors contributing to the 
complexity of shallow water ASW operations is the spatial 
variability of the sea bottom composition. The reflective 
nature of the overlying sediment exerts a significant effect 
on active sonar reverberation levels (RL). Reverberation is 
the sum of all random scattering mechanisms in the ocean 
(ROE C85 2 9b 2237 | In deep water (Figure 1) reverberation 
primarily exists as a combination of two scattering 
mechanisms, volume scattering and sea surface scattering. 
Volume scattering occurs due to scattering from marine life 
within the water column and from inhomogeneities in the 
ocean thermal/density structure. Sea surface scattering is 
due to random reflections (scattering) from a roughened air- 
ocean interface, and is a function of wind speed or wave 
height. 

In shallow water (Figure 2) reverberation from the sea 


bottom typically dominates that due to volume or sea surface 


reverberation. ihe Contribution orebot com sceatteuinamrertae 
total reverberation level (RL) is often so severe that a 
target echo that may have been detected in deep water can 


not be detected in shallow water (Figures 1 and 2). 


D. BOTTOM BACKSCATTERING 

Unlike RL, which is scattering in all azimuthal 
directions, backscattering is only that component of the 
scattered energy that returns to the receiver in monostatic 
(co-located transmitter/receiver) systems. The intensity of 
the backscattered energy from the sea floor that arrives at 


the sonar receiver 1S parameterized as the scattering 


scat 


strength, S, =10 log —2t [Ref. a Dac oo eee |. I is the 
ua 


intensity of the scattered sound per unit area at a unit 
distance and I, is the incident plane wave intensity 
impinging on the scattering area. Measurements of RL and 


See cistian 1G eee ci SCam DCCmEackSCarterca» energy, menave 


indicated that S, is a function of grazing angle (0) and 


bottom type. Backscattering strength traditionally is 
represented using a Lambert’s law relationship where S, = 10 
loge + 0 log sin’ (0). 8 is the grazing angle and Uu 


represents the degree of reflectivity of the sediment 
interface. Based on the deep water data of Mackenzie [Ref. 


Al and. others, sUmick Ref. 3, Boe 2/t-27 3), Gepenks ENnab, ove 


moderately rough topography of varying sedimentary makeup 


and small to intermediate grazing angles Mackenzie’s 


constant (k = 10 log wu) averaged about -29 dB, being 
essentially invariant in frequency and grazing angle. More 
recently, shallow water experiments have indicated an 


assoclaelonmen Mackenzie’ s Constant to the bottom sediment 


type and that the grazing angle dependence may be scaled as 


sin (8) over muddy, absorptive bottoms where volume 


reverberation from the upper layers of the sediment may be 
inporcane wl Ret. = 5] 4 This finding suggests a direct 
correlation exists between backscattering strength and 
sediment type. With this knowledge, it is possible to use 
inverse techniques to infer the nature of the sediment type 


as a function of its reflectivity from RL data. 


E. INVERSION TECHNIQUES 
The "forward problem" in acoustics are those methods in 
which one predicts the characteristics of the acoustic 


propagation field from a knowledge of the governing oceanic 


parameters (e.g., sound speed structure, sea surface 
roughness, geo-acoustic properties, etc.) through which the 
sound is travelling. "Inverse" (or inversion) techniques 


are the exact opposite such that one infers one or more 
oceanic parameters through the analysis of acoustic 


propagation data. [Ref. 6] 


Inverse techniques have a wide variety of applications 
including high resolution geo-acoustic modeling [Ref. 7] and 
Simulation of broadband signal propagation [Ref. 8]. 

Inversion techniques (IT) are used in this study to 
identify features in the beam reverberation time series from 
which bottom sedimentary properties can be deduced. RL 
features such as the slope of the RL time series and higher 
order statistics (HOS) describing the shape of the RL time 


series are analyzed to identify correlations to sediment 


type. 


F. EXERCISE LWAD 99-1 

The Littoral Warfare Advanced Development (LWAD) 
Exercise 99-1 was conducted from 8 to 15 February 1999 in 
the Gulf of Mexico approximately 200 km WNW of Key West, 
Florida, on the southwest Florida continental shelf (Figure 
Bae This test was the ninth in a series of exercises 
performed as part of the LWAD project, a program of shallow 
water exercises designed to identify new technologies which 
have a potential to be transitioned to fleet systems in the 
near future. The LWAD program efforts are concentrated on 
technologies that will aid in the Locate tom and 
classification of diesel-electric submarines and mines in a 
shallow water littoral environment. [Ref. 9] 

LWAD 99-1 included several experiments, one of which 


was for the continued development of inversion techniques. 


The inversion technique experiment had the following 
objectives: (1) record AN/SQS-53C data in a geographically 
variable environment and (2) eollecte ground eeut hh 
information in order to test monostatic and bistatic 
inversion techniques [Ref. 10]. 

The data examined in this thesis was collected during 
event 31 conducted on 9 February 1999 from 0730Z to 1800Z at 
LWAD 99-1 site 1 (Figure 3). The USS Vicksburg ship’s track 
during this event is shown in Figure 4; over 900 pings from 
the AN/SQS-53C tactical- sonar were recorded during the 
event. Event 31 was conducted with a ee present; the 
target track history is also shown in Figure 4. 

le Geo-Acoustics of the South Florida Platform 

The LWAD 99-1 exercise area lies on the southwest 
portion of the continental shelf bordering the west coast of 
Por Gla. The shelf is composed of a karstic Miocene 
platform that is under Pliocene-Pleistocene and Holocene 
sediments [Ref. 11]. The composition of the surface 
sediment was analyzed by Holmes from a variety of sources 
including sediment core samples and seismic data. Figure 5 
1s a depiction of the surface sedimentary distribution at 
LWAD 99-1 site 1 based upon Holmes’ analysis and represents 
the most definitive characterization of the upper ocean 


sediment composition prior to the commencement of LWAD 99-1. 


Howell Hook reef, lying along the western portion of 
Site 1 at a depth of 150 m is a natural border between the 
inner shelf and the shelf slope. Howell Hook is composed of 
coral reef material, while immediately to the east is a 
region of inter-reef sediments. The inner shelf, a region 
of recent sedimentary processes (depositional), contains a 


variety of surficial sediments and is characterized by a 


gentle slope of 0.07°. Farther eastward on the inner shelf 


is a region of foraminiferal sands with sandwaves common. 


Westward of Howell Hook the slope is approximately 0.2° and 


increases through a series of steps until the escarpment is 


reached with slopes as much as 45° [Ref. 12]. 


As part of the environmental characterization for 


Exercise LWAD 99-1, sediment samples were collected. 
However, the results were not available in time for 
inclusion in this study. Preliminary analysis indicates 


that Howell Hook is overlain by a layer of silty sand, while 
the inter-reef region is composed of a clay silt mixture. 
Bast of the inter-reef region the sediment is composed of 
fine to coarse grain sand particles [Ref. 13]. 

Based on the geo-acoustic modeling methods of Hamilton 
and Bachman [Ref. 14], a characterization of the reflective 
properties of the sediments can be accomplished. From Table 


1, the mean grain size can be determined from a description 


of the sediment type. Grain size can then be correlated to 
porosity and saturated bulk density as illustrated in 
Figures 6 and 7. Porosity, which is a measure of the amount 
of voids occupied by water in a material, can then be used 
to determine the compressional wave speed of the sediment 
(Figure 8). The characteristic impedance, the product of 
the sediment density and wave speed, 1S an indicator of the 
reflective nature of a particular sediment composition. 

The geo-acoustic method described above indicates that 
the degree of reflectivity can be determined from the 
sediment type. In particular, reflectivity is proportional 
to grain size. It is anticipated, therefore, that the 
region on the eastern side of site 1 should be the most 
reflective and should have higher reverberation levels. The 
Howell Hook region should be less reflective, but still more 
reflective than the inter-reef region which should be the 
least reflective and have correspondingly lower 
reverberation levels. 

Zs Event Weather and Oceanography 

EO acoustically relevant meteorological and 
oceanographic parameters were essentially invariant during 


the 12-hour duration of event 31, making for near-laboratory 


SOndd talons. The sky was free of clouds with wind speeds 
less than 5 m/s. The sea surface was calm with a 
Significant wave height of 0.5-0.7 m. The surface current 


measured during the event averaged around 0.25 m/s in a 
northerly direction. [Ref. 15] 

A total of 23 CTD observations were taken throughout 
Site 1 between 6 - 9 February. Figure 9 shows the sound 
speed profiles (SSP) calculated from the CTD observations. 


Little spatial or temporal variability between the profiles 


1s observed (< 1 m/s), reflective of the near constant 
oceanic-atmospheric forcing. The mixed layer is seen to 
oscillate between 18 - 25 m. The downward refracting nature 


of the SSP ensures strong bottom interactions of the 
acoustie signal [Ret = 15) 

x Bottom Backscattering Measurements 

Bottom backscattering measurements taken with a 


vertical array during LWAD 99-1 [Ref. 16] indicated that the 


bottom scattering strength (S,) exhibited a sin (@) 


b 
dependence, not the sin’(@) dependence of Lambert’s Law. 


Previous results obtained from a shallow water area south of 
Long Island [Ref. 5] established a sin(@) bottom scattering 
strength relationship with a silty clay bottom having a 
surface sediment compressional wave speed of less than 1600 
m/s. Over sandy bottoms S, suggested a sin’ (9) dependence in 
agreement with Lambert’s Law. 


It is common practice to make backscattering 


measurements in just one location and then extend those 


10 


results to a wider area [Ref. 12], but this practice can be 
misleading. The spatial distribution of sediments in 
shallow water 1s generally highly variable over short 
spatial scale lengths. Scanlon et al. [Ref. 5] showed that 
the spatial scale can grade from hard sand to mud/silt clay 
sediments over distances less than 500 m. 

In this thesis, a bottom sediment map is constructed 
based upon spatially-varying changes in the beam 
reverberation data from an AN/SQS-53C tactical active sonar. 
While the results are considered preliminary due to a lack 
of high resolution bottom samples, they do indicate a high 
degree of spatial variability of sediment types that, only 
on the average, agree with the pre-exercise wide-area 
sediment assessment [Ref. 12]. It will be instructive to 
compare the results of the frequent sediment grab samples 
taken during LWAD 99-1 with the results of this thesis. A 
conclusion that can be stated now is: 

Single point bottom backscattering measurements can not 
predict the high spatial sediment variability and inverse 
techniques must be developed to estimate this variability. 

Another limitation of bottom backscattering analysis is 
the assumption that the scattering occurs at the sediment- 
water interface. Recent studies [Ref. 17] show that the 
scattering occurs from inhomogeneities in the sediment sub- 


bottom for slow (i.e., clay, mud) sediments at low grazing 


—_— 


oe 


angles (Figure 10). Therefore, a full wave model with an 
accurate three-dimensional model of the sub-bottom 
geoacoustic properties must be used to estimate TL in 
shallow water. Just a few years ago it did not seem 
possible to run a full wave model at 3500 Hz, but with 
today’s computing power, it is feasible. Within a few years 
it will be commonplace to run these full wave models in near 
real time and the real time IT sediment-mapping algorithm 


developed here is intended for this future application. 


i 


Bulk 


grain 
Sediment No. Mean grain size Sand, Silt, Clay, density 
type samples mm phi size % % - % g/cm3 
Sand 
Coarse 2 0.5285 0.92 100.0 0.0 0.0 2.710 
Fine 28 0.1638 2.61 92.2 4.1 3.7 2.709 
Very Fine 16 0.0988 3.34 81.0 12.5 6.5 2.680 
Silty Sand 40 0.0529 4.24 57.0 30.9 12.4 2.077 
Sandy Silt 47 0.0340 4.88 30.3 57.8 11.9 2.664 
Silt 19 0.0237 5.40 7.8 80.1 last 2.661 
Sand-silt-clay Zo 0.0177 562 31.7 42.9 25.4 2.689 
Clayey silt 105 0.0071 7.13 7.4 58.3 34.3 2.656 
Silty clay 54 0.0022 8.80 3.9 34.8 61.3 2.715 


Table 1. Continental terrace (shelf and slope) environment; 
Average sediment size analyses and bulk grain densities. 
From Ref. [14]}. 
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II. TREATMENT OF DATA 


A. DESCRIPTION OF LWAD 99-1 REVERBERATION DATA 

The beam reverberation data used in this study was 
collected using the AN/SQS-53C fleet tactical sonar aboard 
the USS Vicksburg, during event 31 of exercise LWAD 99-1 on 
9 February 1999 [Ref. 10]. Beam reverberation data were 
analyzed from two modes of operation of the AN/SQS-53C, the 
SD (surface duct) and VD (variable depression) modes. For 
consistency, all pings that were analyzed used either a 
linear phase modulated pulse or a CP (coded pulse) waveform. 


The SD-mode uses a 300 millisecond pulse length and has 


a 360° azimuthal transmit coverage employing 72 receive 


beams. The VD-mode uses a 500 millisecond pulse length and 
has a 120° azimuthal transmit coverage employing 24 receive 
beams. In both modes, the 3-dB beamwidth is 5 degrees. The 
frequency range is from 2400 to 4500 Hz. The depression 
angle was 3°, while the range scale was set to 20 kyds. 
Event 31 was conducted such that the VD-mode pings were 
directed toward the center of the CG box. This was done to 
get complete coverage of the box and have RL information 


about the same geographic location from different 


directions. When the CG was traversing the western portion 


of the box, on a course of O00°T, the sonar center section 


Ne 


true bearing was 090°T (Figure 11); during the northern leg 
(course 090°T) the sonar center section true bearing was 


180°T (Figure 12), etc. The center section true bearing for 


the VD-mode is the bearing between beams 12 and 13. 


All SD-mode pings (omni-directional) had a center 


Section true bearing of i180°T. The center section true 


bearing for the SD-mode is the bearing between beams 36 and 


8-2 Beam 1 would therefore cover from 000°T to 005°T, beam 
2ELon  O0sen Heomurans etc. Beam 72 would cover bearings 


ot Owes OC ne 


B. HYPOTHESIS 

Reverberation is a combination of backscattered energy 
from within the water column (volume), the sea surface and 
the sea floor. For this study it is assumed that the volume 
and sea surface components of reverberation were constant. 
Therefore, it is hypothesized that for a single ping from 
the AN/SQS-53C different RL beam (bearing) time series 
should only differ because of variations in backscattered 
energy from differing bottom types. This hypothesis is 
further extended to encompass a set of pings (as many as 41 
pings were used in one data set), as long as they occur in 


the same general area and there are no changes to the ping 


type. ; 
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The following assumptions were made to support this 
hypothesis: 
(1) The sound speed profile in the area around the 


source (within ~20 kyds) is constant. 


e As seen in Figure 9, there is little spatial or 

temporal variability of the sound speed profile. 

(2) The output source level (SL) from the AN/SQS-53C in 
any direction does not vary (i.e., the AN/SQS-53C SL is 
independent in azimuth). 

e The difference in performance of the AN/SQS-53C in 

any direction is considered negligible for this 


experiment. Analysis of individual RL time series, 


though, indicates that within 45° of the CG stern an 


increase in ambient (background) noise (AN) due to 
own ship noise is observed (Figures 13 and 14). 
This phenomenon was accounted for in the processing 
of the RL data. 

(3) There are no major topographic changes that would 
significantly change the RL curve (a ne no 
upslope/downslope effects, no pinnacles, etc.) 

e Based on the geologic characterization, there were 

no bottom topography features that would 
Significantly effect the received reverberation 


levels. 
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Ce PROCESSING METHOD 

An example of how the sediment effects reverberation 
can be seen by comparing two different RL time series in 
which distinctly different characteristics can be seen. 


Figure 13 is a single beam (beam 17) RL time series bearing 


082 450a0 One sees that the RL curve reduces to the 


background noise level at around 18 s where the AN level at 


3.5 kHz 1S approximately 35 dB. In contrast, Figure ils 


(beam 37, bearing 182.5°T) shows an example with the same AN 


level, for which the RL curve does not reach the AN level 
Uimtade abou 23.5 Ss. Clearly the beam data is RL-limited 
longer in Figure 15. Tt is suggested that this heightened 
RL is due to the presence of a more reflective bottom. The 
relationship between RL and AN as a function of detection 
range is illustrated in Urick [Ref. 3, Figure 2.2]. 

The general method proposed for RL curve analysis 1S as 
follows: (1) determine a characteristic RL curve for an area 
determined from the average of many pings and (2) compare 
individual beam data to determine its difference from the 
area-wide mean. 

A variety of statistical parameters derived from the 
beam RL time series can be examined to determine if the 
characteristics of the RL curve bear any correlation to 


sediment type. The statistical parameters investigated 
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include the mean slope, standard deviation, skew, and excess 


kurtosis. 
Before eomtinuine, some background about the 
Statistical method used in this thesis is required. The 


Statistics calculated are not based on the entire time 
series but rather from segments of the time series 
partitioned into 0.5, 1, or 2 second increments. In 
addition, a degree of overlap of the data is often useful, 
SUG 2 7m) So mero l@ins 25/50/7755 1Of athe aontervalj2, The 
statistics derived from the overlapped segments are referred 
to as "moving" statistics (moving Se moving skew, 
etc.). For example, a 25-second time series, with 320 data 
points per second, contains 8000 data points. With an 
interval/overlap combination of 1/0.75 the first data point 
of the moving mean would be the mean of points 1 to 320, the 
second point would be the mean of points 81 to 400, the 
third would be the mean of points 161 to 480, etc. Figure 
16 shows an individual RL time series (black) overlaid by 
the 1/0.75 moving mean (white) RL curve. 


Le Example Using Slope Data, SD-Mode, 1/0.75 Moving 
Statistics 


From a data set of 41 pings an area-wide average slope 
is calculated by first taking the 1/0.75 moving mean for all 
of the beams (72) in each transmitted ping. The slope of 


each beam is determined by calculating the gradient (d/dt) 


of the mean RL curve. Calculating the mean from all (72x41) 


no 


Slope data sets generates the area-wide average slope data. 
The area-wide average RL slope for all the pings examined is 
comprised of over 2500 individual RL slope time series. 
Figure 17 shows the area-wide average (solid line) RL slope 
and its standard deviation (--). 

Figure 17 also shows the beam slope for a particular 
beam/ping relative to the area-wide average slope. Each 
point along the single beam slope has an associated range 
and bearing which corresponds to a geographic position. By 
determining the difference between the area-wide average 
slope and the individual beam slope, a geographic contour 
plot can be generated to display the anomaly from the mean 
(Figure 18). Subsequently, different pings can be compared 
to determine if any geographic correlation exists between 
sets of pings to identify a unique but consistent area 
exhibiting similar RL characteristics. 

The other statistical parameters (STD, skew, and excess 
kurtosis) can be analyzed using an identical method to that 


used for the slope. 
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IiIrt. ANALYSIS 


A. ANALYSIS OF RL SLOPE DATA 

ae SD-Mode 

A thorough analysis of all the RL slope plots was 
conducted to create a geographic representation of sediment 
type based on consistent spatial variations in the RL 
slopes. From a continuous sequence of 41 SD-mode pings a 
contour plot was created based on the deviations of the 
individual RL slopes from the area-wide average. 

Through the analysis of the slope data and comparisons 
to that from different pings, RL slope features can be seen 
that are consistent in their geographic location. The RL 
curves from ping 34 exhibit slope characteristics that are 
common in all pings in the set and will be discussed in 
detail. Figure 18 is a geographic contour plot of the 
deviations of each of the 72 azimuthal beams of ping 34 RL 
slope from the area-wide average RL slope. It is 
anticipated that regions where the deviations (from the 
average RL slope) are positive correspond to more reflective 
areas, while regions with a negative deviation are less 
reflective. For example, Figure 19 shows the RL and slope 
data of ping 34, beam 22. At approximately 8.5 s a region 
of increased RL is observed with a corresponding increase in 
RL slope. Figure 20 depicts the region of beam 22 enlarged 


_ 


at 


to show the RL slope details. The same geographic region of 
increased RL can be seen on contour plots of different pings 
(Figures 21 - 23). Note that the RL slope anomaly seen here 
1s consistent in its geographic location andicating that saa 
is due to a geographic feature, and not to a random 
occurrence. 

Figure 24 1s similar to that of Figure 18 but is 
derived from a composite plot of all 41 SD-mode pings. It 
1s the average deviation of RL slope from the area-wide 
average RL slope at each geographic location. Being an 
average of the RL slope deviations, this plot will show RL 
slope features that have geographic consistency, while 


eliminating those features that are not geographically 


consistent. Note that the area within +/- 45° of the CG 


stern were not included in the composite map and accounts 
for the notches, especially to the west, seen in the figure. 

An analysis of Figure 24 shows clear similarities to 
the pre-exercise sedimentary characterization (Figure 5). 


There are 3 distinct geographic locations observed in Figure 


24; east of 84.05°W is the foraminiferal sands area, Howell 


Hook is west of 84.3°W and the middle is the inter-reef area 


which is interspersed with low and high reflective zones. 
The highly reflective zones could be uncovered (or thinly 


covered) reef segments whereas the less reflective spots may 
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be due to depressions/basins filled with unconsolidated 
sediments. It is notable that the foraminiferal sands and 
Howell Hook regions are relatively uniform, while the inter- 
reef region shows high spatial variability of sediment 
types. 

2% VD-Mode 

Due to the limited data available the VD-mode data were 
not used to create a composite map of the exercise area. 
However, individual pings were analyzed to see if there were 
any geographically consistent RL slope features as 
determined from the analysis of the SD-mode data. 

Due to the nature of the exercise the VD-mode pings 
were directed towards the center of the CG box. This was 
done to obtain full coverage of the area within the box but 
unfortunately created limitations for this. study. In 
particular, sedimentary features such as Howell Hook and the 
foraminiferal sands region, which were beyond the boundaries 
of the exercise box, would not visible on the RL plots. 

The VD-mode pings did, however, demonstrate similar 
geographically consistent features as seen in the SD-mode. 
Figures 25 and 26 are individual VD-mode pings from which 
geographic features can be seen. Note that the high and low 
bottom reflection areas observed in Figures 25 and 26 are 
also seen at the same location in the SD-mode composite plot 


(Figure 24). 
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B. HOS ANALYSIS 

The same analysis methods used on the RL slope data 
were employed to identify the correlation of HOS to sediment 
type. Figures 27-29 are plots of RL and the corresponding 
STD, skew, and excess kurtosis. The RL STD and skew plots 
show deviations that can be correlated to RL features. The 
RL kurtosis deviations, however, appeared random and showed 
no apparent correlation to sediment type. The observed RL 
STD and skew features are marginally useful in that they 
only provide an inference that some sedimentary change is 
present. They do not provide eommation about the 
characteristics of that feature (i.e., whether or not it fe 
due to a more or less reflective sediment). 

In order to further investigate RL statistics, a 
distribution analysis was conducted. Figure 30 is a plot of 
skew versus excess kurtosis for the detrended RL beam data 
of “ang 159. The skew and excess kurtosis data shown are 
highly correlated with negative skew (characterized by an 
extended low amplitude tail compared to a Gaussian curve) 
correlated to positive excess kurtosis (characterized by a 
peaked distribution with more data centered around the mean 
compared to a GausSian distribution). Preliminary 
distribution analysis indicates that the reverberation data 
is not only non-Gaussian, but the data spans. several 


distribution types [Ref. 18]. 
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Lewes samportane to note that statisticalseproperties of 
the detrended beam time series are not commonly analyzed in 
classical detection theory applications. Instead, the 
output of a replica correlator is normally the input to a 
Statistical analysis, and normalizations, such as whitened 
Matched filters, etc., attempt to remove the amplitude 
decay, or slope, of the beam reverberation time series from 
the analysis [Ref. 18]. In this thesis we have taken a 
quite different approach by using the beam time series slope 
to invert bottom properties with the HOS of the detrended 
time series as derived product. If the sea objective 
were focused on detection performance, we could have done a 
replica correlation analysis of the residual, or detrended 
beam time series. Instead, we are simply trying to 
correlate the HOS of the detrended time series to the 
environment. Since little is known about the HOS of the 
detrended beam time series, we are just reporting that they 
were asi eh tein non-Guassian and further investigation 
SGnCern1 NomEnelmmEGd Serb tI OnmtineGr1on ls needed) seakrnouch 
the HOS did not correlate well with sediment type, this does 
not imply that they will not correlate with some other 
environmental parameter, such as water depth or sound speed 


DLoLlle . 
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Cc. DISCUSSION 

de Tactical Implications of Sediment Map 

The high variability of bottom sediment type in shallow 
water exerts a strong tactical significance because it can 
unexpectedly affect, often times adversely, the results of 
an acoustic search which is executed under the assumption 
that the bottom sediment is homogenous. Figure 31(A) is 
indicative of a search plan that was executed assuming that 
detection ranges were constant (regardless of sediment 
type). This "mowing the lawn" technique could result in a 
lower probability of detection (PD) than Ba mictnates due to 
regions where detection ranges are degraded because of a 
less reflective bottom. Prior to an acoustic search an 
environmental update could be conducted which includes an RL 
survey from which a sediment map, similar to the one in 
Figure 24, could be used to modify the search plan. 
"Searching by sediments", as seen in Figure 31(B), would 
account for detection range differences due to the bottom 
sediment type. In this case, the acoustic search would be 
executed under the appropriate, spatially-varying PD. 

2% Tactical Significance of Non-Gaussian HOS 

The analysis of the detrended beam reverberation data 
indicated that the backscattered energy is non-Gaussian, 
i.e., the values of skew and excess kurtosis are non-zero, 


and this also has important tactical implications. The 


me 
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Signal processing design incorporated into the AN/SQS-53C 
and other sonar systems is based on the assumption that the 
RL distribution 1s Gaussian. Challenging the Gaussian 
assumption raises two important questions: 

(1) Ts the non-Gaussian nature of the beam 
reverberation HOS due to environmental effects (e.g., 
shallow water acoustic propagation) or because of sonar 
Signal processing design? 

(2) If the RL HOS are due to the environment, will the 
HOS change when calculated for a different environment 
(e.g., deep water, continental slope or soft, muddy 
bottoms) ? 

The first question could be answered by performing the 
same statistical analysis of AN/SQS-53C RL data but taken 
from another environment. If the HOS are similar in varied 
environments, the nature of the HOS is probably due to the 
Signal processing equipment. In this case, the non-Gaussian 
issue would be resolved by making equipment modifications. 

Tf the non-zero HOS are in fact caused by environmental 
factors, then the existing signal processing algorithms, 
from which sonar system design is based, are not valid. 
Currently correlation detectors or .energy threshold 
detectors are optimized, in the probability of 
detection/probability of false alarm (PD/PFA) sense for 


Gaussian beam reverberation statistics. mer the “stonale 4s 
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not Gaussian, then non-Gaussian treatment of the RL data 
must be used to process the reverberation. Solinsky [Refs. 
19-21} has shown that non-zero HOS algorithms can be 
designed which will efficiently mitigate the effects of 
reverberation. 

A fundamental conclusion from this study is that the 
Statistical distribution and its associated properties of 
the reverberation HOS needs to be determined in-situ and 
incorporated into the PD/PFA algorithms to improve echo to 
reverberation processing of the AN/SQS-53C sonar system. 
The HOS characteristics of the environment can be quantified 
through analysis of the abundant AN/SQS-53C data available 
in the Naval Undersea Warfare Center (NUWC) active sonar 
data base. 

The non-Gaussian reverberation statistics problem can 
be addressed in several ways, including the development of a 
neural net (NN) algorithm, schematically illustrated in 
Figure 32. The NN approach can be adapted to the changing 
environment and this technique is highly recommended for 
further development for the AN/SQS-53C. 

Another deterministic approach to addressing non- 
Gaussian statistics is to apply matched field processing 
(MFP) concepts. Here the transmitted signal replica is 
propagated through the shallow water environment uSing an 


accurate propagation model based upon accurate environmental 


ae 


parameters determined by ITs such as the one developed in 
this thesis. Any procedure to process the transmitted 
waveform te accurately reflect the non-Gaussian 
characteristics of the environment will improve’ the 
correlation of the received energy which has been propagated 
in the environment and reflected off. the target and other 
scatterers. If the transmitted waveform is not similar to 
that of the received energy, the correlation will be low and 
detection performance poor. On the other hand, if the 
volume, bottom and surface scattered returns correlate "too 
well" with the transmitted waveform, the active sonar system 
will be overloaded with false targets. MFP or neural net 
algorithms are ways to address the problem of non-Gaussian 
statistics by including the realities of propagation and 
scattering in shallow water in the processing of the 
transmitted waveform prior to correlation with received beam 


data. 
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IV. CONCLUSIONS AND RECOMMENDATIONS 

Exercise LWAD 99-1, which took place in February 1999, 
included an inversion technique experiment during which 
AN/SQS-53C reverberation data was collected. The 
reverberation data was used to investigate the possibility 
of employing inversion techniques to: infer bottom 
sedimentary properties. 

The hypothesis presented in this thesis was based on 
the assumption that the magnitude of backscattered energy 
and its rate of spatial decay were directly related to the 
acoustic reflectivity of the seabed. By examining the slope 
of the RL curves from various sedimentary Reo eeey one 
should be able to identify and associate various reverberant 
areas with unique sediment types. This study examined 
various reverberation level statistics including RL slope, 
standard deviation, skew, and excess kurtosis, to test this 
hypothesis. Of the various statistical parameters the RL 
slope was found to be most clearly aligned with the degree 
of bottom reflectivity. Over less reflective (mud/silt) 
bottoms a decrease in RL slope was noted while acoustic 
interactions with a more reflective bottom (sand) causes an 
increase in RL slope. Based on this finding, a technique 
was developed using the deviation of the RL slope for a 


particular ping/beam from an area-wide average RL slope to 


col 


generate a geographic map illustrating bottom sedimentary 


characteristics. 

Analysis of higher order statistics (HOS), which 
included standard deviation, skew and excess kurtosis, 
revealed no clear correlation to sediment type. Variations 


in standard deviation and skew did provide indications of 
sedimentary change. 

An analysis of the statistical distribution of the HOS 
indicated that the RL data did not derive from a Gaussian 
population. It also showed that the distribution did not 
appear to come from any well-known distribution (e.g., 
Poisson, log normal, Rayleigh, etc.). The non-Gaussian 
nature of the RL data has serious implications for sonar 
Signal processing because most Signal processing algorithms 
assume a Gaussian distribution to the noise field. 

It is recommended that the RL slope technique for 
generating a sediment map be used by fleet assets as part of 
efforts to develop in-situ and archived bottom reflection 
databases and that these results be incorporated into 
shallow water ASW search planning. iy “adcaitaeom,  fulet hers 
investigation into the nature of RL statistics and its 


impact on sonar system performance should be undertaken. 
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APPENDIX A. FIGURES 







Combination of surface and 


/ volume reverberation 


Target Echo 


Time/Range 


surface scattering 





Volume Reverberation D ae . 
AY ater 


Figure 1. Schematic and characteristic RL curve for a 
deep water environment. Reverberation is due primarily 
to backscattering from volume and sea surface scatterers. 
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Combination of surface, volume 
Pa and bottom reverberation 






Time/Range 


surface scattering 





ANAANANANSN NS 
Shallow Water 


Figure 2. Schematic and characteristic RL curve for a 
shallow water environment. Reverberation is generally 
louder being a function of volume, sea surface, and 
bottom backscatterers. High RL levels often mask the 
target echo. 
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Figure 4. Event 31 experiment geometry. Over 900 pings 
of the AN/SQS-53C tactical sonar were completed between 
0730Z and 1800Z on 9 February 1999. From Ref. [10]. 
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Figure 6. Mean grain size versus porosity. Round dots 
represent continental terrace (shelf and slope) samples. 


From Ref. [14]. 
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Figure 7. Mean grain size versus saturated bulk density. 
Round dots represent continental terrace (shelf and 
slope) samples. From Ref. [14]. 
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Figure 8. Porosity versus sound velocity, continental 
terrace (shelf and slope). From Ref. [14]. 
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SITE 1 CTD DOWNCASTS 
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Figure 9. 23 total CTD downcasts were taken at site 1 
from 2/6/99 to 2/9/99. Bottom interaction of acoustic 
energy is ensured due to downward refracting nature of 
sound speed profile (SSP). From Ref. [15]. 
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Figure 10. Scattering mechanisms are from (1) the water- 
sediment interface, (2) sub-bottom volume, and (3) the 
sediment-basement interface. From Ref. [17] 
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Figure 12. VD-mode ping with center section true bearing 
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the CG box. 


Ping: 189, Mode: SD, Beam: 17, Brng{T): 62.5 
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Figure 13. RL time series for a forward projected beam. 
RL reduces to the ambient noise (AN) level of 
approximately 35 dB at about 18s. Compare with Figure 14 
(stern beam) to see difference in AN level. 
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Ping: 189, Mode: SD, Beam: 54, Brng(T): 267.5 
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Figure 14. RL time series of stern beam. The ambient 
noise (AN) level is approximately 20 dB greater than non- 
stern beam due to ship noise. 
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Ping: 189, Mode: SD, Beam: 37, Brng({T): 182.5 
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Figure 15. RL curve reduces to background noise level at 
around 23.5 s. Compared to figure 13 the longer RL decay 
time is associated with backscattering from a more 
reflective bottom. 
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Ping: 189, Mode: SD, Beam: 17, Brng(T): 82.5, Interval: 1, Overlap: 0.75 
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Ping: 189, Mode: SD, Beam: 17, Brng(T): 82.5, Interval: 1, Overlap: 0.75 
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Figure 17. Slope of the 1/0.75 moving mean (-*) for an 


individual beam compared with the area-wide average slope 
(solid line) and its standard deviation (--). 
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Ping:0034 Mode: SD, Interval: 1, Overlap: 0.75 


CS 
—_ 


' 
' 
4 
4 
‘ 
4 
4 
‘ 
‘ 
4 
i 
’ 


2a. 


(cass) adojs Ty ebevare Woy uoieAeq 


Jano 

Zo. 
25.65 
20105 


(N) apne 


. @ 
dl 
a we L 
ry 4 
s 4 


-10 


63:95 83.9 


B4 


meeedtowe ween dane caceae 





84.05 


Longitude (VV) 


84.3 84.25 842 8415 84.1 


Geographic contour plot of 


Ping 34. 


Figure 18. 
deviations of all 72 beams of slope data from the area- 
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reflective sedimentary areas than the cooler colors. 


wide average slope. 
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Ping: 034, Mode: SD, Beam: 22, Brng({T): 107.5, Interval: 1, Overlap: 0.75 
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An increase in RL at 8.5 s 


beam 22. 


corresponds with an increase in RL slope. 


Ping 34, 


Bagure 19. 
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Ping:0034 Mode: SD, Interval: 1, Overlap: 0.75 
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Figure 20. Ping 34, beam 22. Enlarged contour plot to 
display the region of increased RL. 
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Figure 21. 


geographic location as Ping 34. 
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Figure 22. 


geographic location as Ping 34. 
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Ping: 35 | Mode: SD, Interval: 1, Overlap: 0.75 
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LWAD 99-1 Sediment Map 
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Figure 24. Composite plot of deviation of RL slope from 
the area-wide average for 41 SD mode pings. Warmer 
colors represent more reflective sedimentary areas than 
the cooler colors. 
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Figure 25. VD-mode ping 2. Note that the features 
pointed out are also seen in the SD-mode composite plot. 


of 


5 
25,7 sli , — 


—2 


25.7 | 
z= 
@® 
3 10 
: 
25.65 |-------- 


Deviation fram average RL slope (S92) 


SS ——— 
- 2a [ 
| ec ~ h 
f~ . NA 
A i NT 
f — 
f \ A 
Jy \ 
, “Y ~~ » 5 L 
=o - 


256 





Longitude (VV) 


Figure 26. VD-mode ping 4. Note that the features 
pointed out are also seen in the SD-mode composite plot. 
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Ping: 034, Mode: SD, Beam: 23, Brng(T): 292.5, Interval: 1, Overlap: 0.75 
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Figure 27. 
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Ping: 034, Mode: SD, Beam: 23, Brng(T): 292.5, Interval: 1, Overlap: 0.75 
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Ping 34, 


Figure 28. 
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Ping: 034, Mode: SD, Beam: 23, Brng(T): 292.5, Interval: 1, Overlap: 0.75 
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"Mowing the Lawn" “Searching by sediments" 


Figure 31. Old search method (A) vs. new search method 
(B) with knowledge of local sediment types. 
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Figure 32. MTwo-Dimension feature space. From Ref. [21]. 
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APPENDIX B. MATLAB FILES 


The matlab program files enclosed are for analyzing 
beamform data from exercise LWAD 99-1. The program files do 
not include the actual beamform ping data, or the ping 
parameter files. The matlab files included here are for 
guidance in creating a graphical user interface (GUI) for 
analyzing beamform data. 

Beamform data files contain the following variables: 
bandwidth - bandwidth (Hz) 
bm£ - complex beamform data 
centerfreq - center frequency (Hz) 
ping_date - date (e.g. ‘'02/09.99') 
ping _time - time (e.g. ‘11:18:15’) 
pulselength - pulse length (ms) 
samplerate - samplerate (Hz) 

The parameter file ‘evt31l.mat’ contains a vector of 
times and corresponding CG and TGT positions from LWAD 99-1 
event 31: 
evt31_ time - time in seconds from midnight 
evt31 date - event 31 date (e.g. ‘02/09.99’) 
cg_lat - vector of CG latitudes 
cg_lon - vector of CG longitudes 
tgt_lat - vector of TGT latitudes 


tgt_lon - vector of TGT longitudes 
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The parameter file ‘’evt3lprm.mat’ contains a vector of 
times and corresponding ping parameters: 


botdepth - bottom depth (fm) 
course - CG course (°T) 
ctrbrng - Center of ping section bearing (°T) 


deangle - sonar depression angle 
speed - CG speed (kts) 


timeday - time in seconds from midnight 
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dP dP dO dP dP dP oP AP DO dP dO AP DO oP DO DO dP dP AP GP oP oP DO dP dO DP dO DP dP dP GP dO dP GO oO cP oP dO dP GO oO oO oP oP oO oP 


Run a 


GUI to analyze LWAD 99-1 data 


A. Begin by selecting a ping to analyze 
The following will occur: 


aly 


A map of the LWAD 99-1 exercise will be created 


showing the following: 
(a) blue cg track history, red circle is cg 


current location 


(bys redmtqtretrack Histomy,, blue circle is tgt 


Current elecat lon 


(c) ping section is shown in yellow with dashed 


yellow center 
section true bearing line. Yellow asterisks 
show 10nm range. 


(d) latitude/longitude grid lines (black dashed 


lines) 


(e) Sediment map showing Howell hook (black circles) 


(2) 


oF 


with phi size approx 4-7, clay silt (black dots) 
with phi size approx 6-8, and sand 

(black slashed lines) with phi size approx 1. 
This map is based on data received from 

J. Fulford (NRL Stennis). 


A list of ping parameters will be displayed in 
the ‘LWAD Parameters’ section. 


A data analysis section will be created 
in the GUI. 


create a contour plot press ‘Contour Plot’ button. 
This plot is created by determining the 
difference between area averaged slope data and 
slope of individual beam data. 

The plot is a filled contour plot. 

A regular(non-filled) contour plot may be 

created by using the function ‘lwad_contour.m’. 


Note: when the ‘contour plot’ button is pressed, the 


contour plot is created 
uSing an interval of 1 second and an overlap of 


.75. The interval/overlap is not what is selected 


in the LWAD Data Processing Program GUI. If you 
want to investigate different interval/overlap 
combinations, you’ll probably want to run 


‘lwad_contourf’ function separate from the LWAD 
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dP dP dP dP dP dP dP DP oP ADP DP dP AP dP DOP dP AP dP dP dO dO dP dO OP dO GP GO dO dP GO ae dO dO dO dP oe? 


dP dP dO dO oP 


Data Processing Program GUI. 


C. To analyze beamform data 


(1) select the type of plot(s) to display 
- 20log|data| - is RL data in dB with averaged 

data in red 

- Slope - slope of RL data in blue with area 
averaged slope in red. Dashed red 
line is +/- 1 standard dev. 

- STD - standard deviation in blue with area 
averaged std in red. Dashed red line 
1s +/- 1 standard dev. 

- Skew - Skewness in blue with area averaged 
skewness in red. Dashed red line is 
+/- 1 standard dev. 

- Kurtosis - Excess Kurtosis in blue with area 
averaged excess kurtosis in red. 
Dashed line is +/- 1 std. 

(2) Select the beam you would like to analyze 

1-24(VD) or 1-72(SD) 

(3) Select the interval and overlap to be used by 
mmean.m, mstd.m, mskewness.m, mkurtosis.m. I 
usually use 1 sec and .75 overlap. Other 
combinations may be used to see the effect. 

(4) Press ’Finish’ to see the plot(s). 


(5) Press ‘Close plots’ to close all plots but the 


event map, and the LWAD GUI. 
D. To analyze another file press ‘Browse’. 


E. Press ‘Quit’ to end your session and close all 
open plots. 


See also LWAD_CONTOUR, LWAD_CONTOURF, LWAD_MAP, 


LWAD_ POSIT, LWAD_PRM, MMEAN, MSTD, MKURTOSIS, 


MSKEWNESS. 


Created: May 8, 1999 
David Schalm 


bast Moditived-: July 26, 1999 
David Schalm 


menu_num = 1; 

figure 

load ‘e:\schalm\lwad 99-1\data_drive’ 
lwad_menul 
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S$LWAD_MENU1 


% 

% THiS Wiel "ereste a Ul to select a £1le for the LWAD 
% program 

% 

% See also LWAD 

% Created: April 29, 1999 

% David Schalm 

% 

% Last Modified: July 26, 1999 
% David Schalm 

=e weicu. 

Set (fl1g es 


‘Wits = Doints >. 

 Groplenad (0. S08 8 Oe 8], 

‘MenuBar’, ‘none’ ans 

‘NumberTitle’, ‘off’, 

‘Name’, ’LWAD Data Processing Program’, 
‘Position’, [450 240 300 ee 

‘Resize’, POLE! i 

HECOntco mC Unaks! points” : 

PERC LG OMe ISL PAOine 0. 8 on B14 
JHotTtzZontalAlionment 5’ Lert’ | 
SPOchedlom = h0. 260 180-10 ]., 
‘String’, ’LWAD Data file’ 
‘FontWeight’, ‘bold’, 

1 SEVLe  ECXE Ih; 

O° — Ureonmwrol Unies | DOolnts ” 
‘Backorounacolonm 7 ll il], 
‘HoruzontatlAli onment’ ;’ lett, 
‘Positron -@rlo 260) 160 15), 
‘Stung ,~,tdatvaLaraye "* mat’ |i, 
Seve" edit" )- 

callbacki = [’ [data_file,data_drive] = ' 


LUlGeemrle (data ldrive, (6x kets ae 
‘set(fig,’ ‘’’Name’’’ ‘“,’ on sci Pacey, Nee ee 
De Tee Gatamedle ==.0, load : 

ae eee ae 99-1\data_drive’’’ ‘';"., 

‘set (fig, Name’ 90", * 

‘''TWAD Data Processing ama aoa 


‘,else, lwad_menu2, end’]; 
Wiconerol(: Unies = ports’ 
‘Callback’,callbackl, ... 
YPOsteLon sale 00260540 15), 
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‘String’ )’ Browse ae 
callback2 = [’ButtonName=questdlg(’ 
‘’*’Are you sure you want to quit?’’’... 
1 MEMEO Ee SEA a ON ee 
1 EON OC a ee BS ee ee 
‘if length(ButtonName) == 3, close all, end’); 
Wacontrol( /Undteewe De isa isu 
‘Callback’,callback2, ... 
“POSitiom oil. 0 0 aes |. 
(CS tere Gh oe POM Minte Nae 
LTLeConerol sUntES + pedis. 7 ea 
“Backorounacolor fiOs5 02580751 
‘HorizontalAlignment’,’left’, 
‘Position’, [5 250 290 1], 
COErING 3 0 aon 
‘Style’,’text’); 


ae) 


LWAD_MENU2 


This will display the file parameters, and add a UI 
section for selecting plot(s) and plot options. 


Note: when the ‘contour plot’ button is pressed, the 
contour plot is created using an interval of 1 
second and an overlap of .75. The interval/overlap 
1s not what is selected in the LWAD Data Processing 
Program GUI. If you want to investigate different 
interval/overlap combinations, see 
slwaGuGOueoUure-l- 4s prOU li prebably want..co run 
‘lwad_contourf’ separate from the LWAD GUI. 


0 dP dP dP dP dP BP cP AP DOP oP oO oP oP GP? 


See also LWAD, LWAD_CONTOURF 


% Created: April 29, 1999 

% David Schalm 

% 

% Last Modified: July 26, 1999 
% David Schalm 

fac. = Oem, 


save ‘'e:\schalm\lwad 99-1\data_drive.mat’ data_drive 

fname = [data_drive data_file] ; 

load(fname, ’-mat’ ) 

[course, speed, ctrTbrng,botdepth,deangle] = lwad_prm(fname) ; 
fFeglat, colen, tgqtlat, tgtlon) = lwadypesit (iname) ; 
num of beams = size(bmf,2); 

set (h0, 'string’ ,data_file) ; 


lil CON EBON Unites 2 potmcG = 2. 
Boer GcomaCo lon PLOws 0. 8 0. ly 
On one camer cere 
POsmeron 1.0250 1e0e10 |e 
‘String’, ’LWAD Parameters’ 
‘Fontweggne., “bolas, 

‘Style’, ’text’); 

LEONEL OM Umites Dol: S 
—Peaer aroun color eo 6 0. 8 0. 8], 
‘HorizontalAlignment’, ‘left’, 
‘Pesteaon sto 20 140 10)7 
‘String’, [’Ping date: ’ ping_date], 
‘Stvle’> “texe): 

HLeEGntero (“Unies oe poles’ |  . 
Pomerarcunceolon ROG 0. 8 0. 8], 


ial 


HOrrZonta lAigmmecmits bert 4 
“POSItiOn ol LS 0 e200 si 4Oa 0); 
‘String’ jl” Ping sieme- 6. ep ping time, 
‘'Stvle" > text Ve 


lwad_map (fname) ; 

figure(rad}) 

UlLGOnErOl( Unitus, polnes : 
pee he rounaeolont =O. Ss 0. 8 0. Sule 
HOE zZOnkalldaemement |, lett > 
“POSiETON , (LO Z2ZOCNI40 1.0) , 

“String MCG Posteion: 2 mum2s cr (ecqlatyme ine 
numZstr@eoqglon) Ww )-; 
(Style =. Gext > 

U1 CGOmemem Units’; points. ee 
Pee ounce lor Os 0. 8 0. Sly 
‘HorizontalAlignment’,’left’, 
"Poscttrom = | 50) 200 140 10]; 
‘String | i PGP Posieion= “ mumZ2cer (tctibom jell 

nMuNZSES( EGEen) “Wo |e 
SEV Le “ECE. 

HlcomMenol, Wntes ). Doles. ne 
WBackoBelmMacGlLor - | 0.6 sm em0 . 6.1. 
‘HorizontalAlignment’,’left’, 
*POStELON: moe 90m 4 0 Olea. 
‘String’,[’CG Course(T)/Speed(kts): ' 

num2str(course) ‘/’ num2str(speed)], 
"Sty Le" gos Geer 4). 

ULCOnEroO!] ( Unite ) pomnes. : 

Backer oundeolon: Ose 0. 8 0. Sa 
‘HorizontalAlignment’,’left’, 

POSPeEMom ells 0190; 140. 10), 

String nl Cce@eecre Searing: “ snumZzstr (Gtmlornc) iy, 
‘Style’, ’text’); 

WECONETOIN Undis 3) polnts |e 
acer suneolor. AOE 0. 8 0. 8], 
‘HorizontalAlignment’,’left’, 

PPOSPTETON Sal Ol SO a0 AO, 


“String’ , | Bottom Depth(im): ~ numZstr (boredeorn) |; 
‘Style’, text:); 
Wileotero.l | Units 3 pores: 


Pew OE nacolloe OS 0. 8 0. Sule 
‘HorizontalAliognument,,’ lete 
‘Position’ ,;,[(150 130 =eaCeao, 
‘String’, [’Depression Angle: ’ num2str(deangle)], 
Sty le text) 

UicOnerol Units | pedmrc. 


12 


“Bachemoutnicwenkem. (0.8 0.85 0.8); 


‘HOorazeonralealigument’, left’, 
‘Position vam t70 140 10), 
‘String’, [’Bandwidth: ‘' num2str(bandwidth) ’ Hz’], 
Sty lee Abext ):-> 

ULCOMNEEOINat es, Points” : 
 BackoGwoumnaGovor , (02507 c 0.6 |) 
‘HOrizZzonealAlagnment~, left’, 
POsSi tiem er oo sl 70 140710] 
‘String’, [’Center Freq: ‘' num2str(centerfreq) 

Peiaeleocy le’. text 4); 
Ulakexehelcieto IU) ° ha kets Ponca: 


Pee olor! alr S 0. 8 0. Sale 
‘HorizontalAlignment’,’left’, 
UPOsieeadton i bOmmieoO 1 14G> 10], 
‘String’, [’Pulse length: '’ num2str(pulselength) 
menses Style ) text”): 
MECOnE he 4 4Umitswe. POlmts. oman 
Peace anotndcolon FiO 26 0. 8 0. SB], 


Pe ence onnene. lanes 

Looavtron: sil b5s0 160° 140 10), 

‘String’, [’Sample rate: ’ num2str(samplerate) 
vpuizeale. Style. “text ) 


% Create time vector 
cleies, Sigenee 2 ee 


data_duration = length(bmf)/samplerate - (1/samplerate)... 


+ data_start; 
time = (data_start:1/samplerate:data_duration) '; 


WICOntrel (aun Es spelnmtes: : 
Backaroundeolor? Oke ts 0. 8 0. cue 
‘HorizontalAlignment’, ‘left’, 
‘Poss taeone., He aS 40 LO |; 
‘Steing - (| Data, Duration: " num2sitr (datamciicat rom)... 4 

‘ Sec'], 
‘Style {7 texte Ms 


Va COnLrOl sUmilis wa eolnts 
Backhoe enon: PleOu. oS 0. > 0. alr 
‘HorizontalAlignment’,’left’, 
ae Ti Seo 290). 


“String, 
Stevie” Bates : 
iconenone Units. points’ 


Bee oie onlora mLOSS 0. 8 0. Slr, 
Hora zZentaltAlaonment’ ,’ lekt’, 


1S 


“Pesteron’ ,.). 0 ie Ca OO 


Toe Ge : 

‘select the plots ® you would) fikew permits eiay 

‘FonkEwelgne. =. bola, 

“Style CeExt a) 
hie—-uaigentirsete Units a points” : 

PBaG kgEounaColer Salmo O60. Se 

‘Hori zontalAlagnmentr’ > lett 


“PoOsme Lom caleo OO moo.) 10 |; 
“String mecO log cata, ’ 
‘Style’ “eneckbox! ):- 

h2 = Geeeneror DnTes ] points * dae 
“Backdrounadcolon - tO. 8 0.8 0. 8], 
‘HorizontalAlignment’,’left’, 
“PosteLon sallOrs > £O0 LO, 

(SELING aS lope” . 
‘Style’, pone i" 

ho = swtconeroliwwmbrs.) polmra: Ma 
 Backemouncdc Omow s LOse O25 JO026 |) 
‘HorizgontalAlignment’, lett’, 
“Positron welo 70 60 0], 

Ovengaoregae, | SUUB 
‘Style’, ‘checkbox’ ); 

ha = mconcrom et masks “3 Dodmnes * Oia 
“Backo LOunacoloer - Soe 8 0.8 0. Br 
‘HorizontalAlignment’,’left’, 
/Poscwimon’ (0° 20-60-10 |x 

"SEring i wnt oslS ” 
‘Style’, ‘checkbox’ ); 

NS- = ONConkuema Umirtes.) points” as 
‘BackareunaGolor’ 10.8 0.8 0.3), 
‘HorizontalAlignment’,’left’, 
‘Post elton O55 60-201, 

‘String’, ’Skew’, : 
‘Style’ ,eechneckbox’)% 

ULCOnErol( Unmmes | polnes 
‘BackGrounadcolee e(0.6 C2e80- 5), 
‘HorizontalAlignment’,’left’, 
oct iton’ , 65 L00 (65 210 |., 
‘String’, ’Beam Selection: ’ 
‘Style’; Cex |); 

beam = {1:1:num_of_beams}; 

hn6 =sUtcoenerol (7Uinikes 7 DOImntes. 
‘HorizontalAlignment’,’left’, 
“BackgroundColer | hl. lana 
“P@s@eton ess 420iGe. 55 ie 
‘String’ ,beam, 
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‘Sty lever i 1S toes 
Atagoic ston dveje4 Jj 
UUECONTCEO Men ris. Dolnts' , 
i paerorounacodor: PaO 0. 8 0. 8], 
Mori cone anicnment® lett’. 
(POSitiememeneo®. 100565 10), 
‘Strings. Interval (sec)’, 
“Sty \e a eeeexk) 
Pieter valmse besa. a) LZ) 
hy =] Uateomane hi Unies’. pormes - 
‘Backgmoumacolor’ 7(1 1-1}, 
‘HorizontalAlignment’, ‘left’, 
POSiteronm wll6O 40660) 5515 
String eunterval sel, 
LSEyre . “1istbpes 
WULCoOnissOl mit S 6.) Dowmes 7 
Peccgremndcelor Oem Oe-s 0. ae 
‘HorizontalAlignment’,’left’, 
APeSite Lom .(49 0. 200: 65,10), 
CS Meer Over | apeO— 1)", 
“Style we. EEXt” )- 

Over lapse! — 025) 35 75): 
how UICoOmeno FOsumines “Dornts ’ 
(BeckoQmouna@ovor sab ie 17 
‘HorizontalAlignment’, ‘left’, 

~POSTEVoOn sacs 0.40.60 55], 
‘String’,overlap_sel, 


‘Style’, ‘listbox’'); 
Calleack4 = [set (iig, Ss Name’ mee . 2. 
a +GeEtCIingm@ontour Plotm..'%" *) >’ 
elwad Conteuret(iname, 1, .75) >’ 
‘set (fig, ’ ‘i’ 'Name’’* aera 
‘''TWAD Data Processing eeoorame aa 


ULCOnErcIN Units — poants 
‘Callback. See 
CPosth=tom eel SO 10 60) LS): 
Cling ee Contour lot, jr 

callback3 = [’set(fig,’ ‘'’Name’’’ ‘’,’ 

PET Ceteing. Plots. oC y=? 
‘lwad_menu3’]; 

HrCOntL LOL (Unt eSs ~! polnes” 
‘Callback’, callback3, 
POs Ledon ells Ohm 0 15), 
COE VAM) 4) Beene in 9) 

set(fig, ‘Name’, ’LWAD Data Processing Program’ ) 
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$6LWAD_MENU3 

% 

% This script file is used to create the selected plots. 
% 

% See also LWAD 

% Created: April 29, 1999 

% David Schalm 

% Last Modified: July 26, 1999 

% David Schalm 


beam_selection = beam{1}(get(h6,’value’)); 
datamtiriler—eget (nO, string’ )- 


interval = interval_sel{get(h7, ‘value’)}; 
overlap = overlap_sel{get(h8, ‘value’) }; 
spacing = interval - interval*overlap; 


N = round(interval*samplerate) ; 


% Create time vector 
[mtime,ao] = mmean(time,N,overlap) ; 
1£ ao ~= overlap 
overlap = ao; 
disp([’Note: actual overlap is ’ num2str(ao)]); 
end 


666% PLOT fue lene DATA SESEESEESEESES 


S% 
Display plots of data/mean, slope, kurtosis, skewness 


oP dO 


if datalftile(s 76) e== “Sp: 
area_stats_fname = [’e:\schalm\lwad 99-1\area_stats\’... 
data_file(5:6) ‘_stats’ 
num2str(interval) ’.mat’]; 
else 
area_stats_fname = [’e:\schalm\lwad 99-1\area_stats\’... 
data_file(19:20) ‘_stats’ 
num2str(interval) ‘.mat’]; 
end 


load(area_stats_fname) 


if get(h8,’‘value’) == 


area_stats = stats_0O; 
elseif get(h8,’value’) == 2 

area_stats = stats_25; 
elseif get(h8,’value’) == 3 

area_stats = stats_5; 


a6 


elseif get(h8,’value’) == 4 
area_stats = stats_75; 
end 


log_data = 20*1log10 (abs (bmf(:,beam_selection) )); 
thetal = ctrTbrng - num_of_beams*5/2 + 5/2; 
theta = thetal + 5*(beam_selection - 1); 

if theta < 0, theta = theta + 360; end 


ticle Stwtncme sl Fing-s soata f£1le(l4:16)m", Mode: SD’ 
‘, Beam: ’ num2str(beam_selection)... 
 eBEngtt).. 7 Num2str (theta 
‘, Interval: ’ num2str(interval) ’, Overlap: ' 


num2str(overlap) J; 


% plot 20log|data| 
if Cet (mie volte. )i.——o 
ENGUae woe — OCT ; 
plot(time,log_data,’k’), title(title_string), 
ylabel(’Reverberation Level (dB)’), xlabel(’Time (sec) ’) 
axis POmeeta@elrabione. | 50) eegridrons. nodcar om 
[mean_data] = mmean(log_data,N, overlap) ; 
plot (mtime,mean_data, ‘w’ ) 
end 


% plot the slope of 20log|data| 


1f get(h2,’‘value’) == 1 
figure soma ek - 
[mean_data] = mmean(log_data,N, overlap) ; 
slope_data = gradient (mean_data) ; 
plot (mtime, slope_data/spacing, ‘k*-’),title(title_string) 


ylabel(’RL Slope’),xlabel(’Time (sec)’),grid on, hold on 
plot(mtime, area_stats(1,:)/spacing,’k’) 
plot(mtime, (area_stats(1,:)-area_stats(2,:))/... 


spacing, ’--k’) 
plot(mtime, (area_stats(1,:)+area_stats(2,:))/... 
spacing, ’--k’) 


legend([‘’Beam ’ num2str(beam_selection)... 
‘ Slope’], ‘Average Slope’, ‘Average Slope STD’) 
end 


% detrend the data w/ a 3rd order polynomial 
log_data = log_data-polyval (polyfit(time,log_data,3),time) ; 


% Plot the STD 


it Cer (ios ava like). —— ai) 
Figure. 9S) = OCr. 


Cu 


[std_data] = mstd(log_data,N,overlap) ; 

plot (mtime,std_data,’k*-’), title(title_string) 
ylabel(’STD’), xlabel(’Time (sec) ’) 

Gridy,on, hold won 


plot(mtime, area_stats(3,:),’k’) 
plot(mtime, (area_stats(3,:)-area_stats(4,:)),’--k’) 
plot(mtime, (area_stats(3,:)+area_stats(4,:)),’--k’) 
legend([{[’Beam ’ num2str(beam_selection) ’ STD’],... 
‘Average STD’, ‘STD of Average STD’) 

end 

% Plot the Kurtosis 

Lf get (ha4nevalue’) == 1 
EPG e «sou —FOCt 
[kurtosis_data] = mkurtosis(log_data,1,N,overlap) ; 
plot (mtime, kurtosis data, ’k*-’),title(title_string) 
ylabel(’Kurtosis’), xlabel(’Time (sec) ’) 
GLa on, nold om 
plot(mtime, area_stats(7,:),’k’) 
plot(mtime, (area_stats(7,:)-area_stats(8,:)),‘’--k’) 
plot(mtime, (area_stats(7,:)+area_stats(8,:)),‘--k’) 
legend([’Beam ’ num2str(beam_selection) ‘’ Kurtosis’],.. 

‘Average Kurtosis’,’STD of Average Kurtosis’ ) 

end 

% Plot the Skew 

if get(h5,’‘value’) == 1 
PUGur ere) = “Ger: 
[skew_data] = mskewness(log_data,N,overlap) ; 
plot (mtime, skew_data, ’k*-’),title(title_string) 
ylabel(’Skew’),xlabel(’Time (sec)’) 
dria on, Nola Jon 
plot(mtime, area_stats(5,:),'’k’) 
plot(mtime, (area_stats(5,:)-area_stats(6,:)),’--k’) 
plot(mtime, (area_stats(5,:)+area_stats(6,:)),’--k’) 
legend([’Beam ’ num2str(beam_selection) ‘’ Skew’],... 


‘Average Skew’, ’STD of Average Skew’ ) 
end 


set(fig, ’Name’,’LWAD Data Processing Program’ ) 
figure (fig) 
Wee niteaote(” Unises. se se @en ates sets 

ZC alilbackuel closeout frac  £0O))) "ih 

BE oswtdeil 7 (65 whOmo eles |), 

YGering ,’ Close: Plees-)-: 
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function [area_stats] = lwad_stats(file_type,interval... 


,overlap) 


‘SLWAD_STATS Determine the statistics for the LWAD 99-1 


dP dD CO oP BG CP ADP CP DO DO CO CO CO cP DP DP DP DO CO CO DO DO DP CO CP GH GO DO dO DO oO dC GC HO CO cP dO CP CO oP oP oO oO GC? 


exercise area. 

[AREA _ STATS, ERROR_LOG] = LWAD_STATS (FILE_TYPE, 
INTERVAL, OVERLAP) will return the statistics for the 
selected parameters. The function will calculate 
Statistics for each individual beam, average these 
statistics, and determine the standard deviation. 


FILE TYPE aS either ‘sd’, ‘cp’, or ‘cw’ ... depending 
on which mode you want to analyze. 


INTERVAL and OVERLAP are required by funtions MMEAN, 
MSKEWNESS, and MKURTOSIS to calculate the respective 
Statistics. If no INTERVAL or OVERLAP are specified, 
the program defaults to INTERVAL = 1, OVERLAP = .75. 


The mean and std for the following statistics are 
determined: 

(1) Slope 

(2) STD 

(3) Skewness 

(4) Excess Kurtosis 


Note: STD, Skewness, and Excess Kurtosis are determined 
from ‘detrended’ data. The data is detrended by 
subtracting a best fit 3rd order polynomial from 
the data. Slope is not detrended. 


The output file AREA_ STATS 1S organized as such: 


Row 1: mean of the slope 
Row 2: std of the slope 
Row 3: mean of the std 

Row 4: std of the std 

Row 5: mean of the skewness 
Row 6: std of the skewness 
Row 7: mean of the kurtosis 
Row 8: std of the kurtosis 


The number of columns will be a function of the 
selected INTERVAL and OVERLAP. 


If the time series data contains a zero value this will 
cause errors in the statistics because log(0) = -inf. 
All beams which contain a zero value are discarded. 

An error log may be obtained as an output. 


eo 


Created: May 25, 1999 
David Schalm 


Last Modified: May 25, 1999 
David Schalm 


oP db dO dO oP 


load data_drive 
drive_num = data_drive(l1); 


ieee ehior = — 0 

file_type = input(’Select file type (sd,cp,cw)...’); 
end 
if G2lemevyoe == sd’ 


filepath = {[drive_num ... 
‘:\LWAD 99-1 Data\Tape 17\sd\‘] ... 
[drive_num ‘:\LWAD 99-1 Data\Tape tonsil” ee 
elseif file_type == ‘cp’ 
filepath = {[drive_num ... 
‘;\LWAD 99-1 Data\Tape 17\vd\cp\'] ... 
[drive_num ‘:\LWAD 99-1 Data\Tape 18\vd\cp\']}; 
elseif file_type == ‘cw’ 
filepath = {[drive_num ... 
‘:\LWAD 99-1 Data\Tape 17\vd\cw\'] ... 
[drive_num ‘:\LWAD 99-1 Data\Tape 18\vd\cw\']}; 
end 


row_index = 0; 
error index = 0; 
Lor yb: =" 

clear fname 

D = dir(char(filepath(b))); 


[fname{1l:length(D)}] = deal(D.name) ; 
%* get rid of ’.’ and ’..’ filenames 
fname = fname(3:length(fname) ); 
Lox a = 1:1:length(fname) 
oad([char(filepath(b)) char(fname(m) ) ]) 
sie read file is ’ filepath{b} fname{m}] ) 
x =" (1: length tbomt)) 
N = round(samplerate*interval) ; 


HUM Onebeams = Size lomi,2) 


% Create ecime vector 
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data_start = 1; 
data_duration = length(bmf)/samplerate - 
(1/samplerate)... 
+ data_start; 
time = (data_start:1/samplerate:data_duration) '; 
[mtime,ao] = mmean(time,N,overlap) ; 


for ny =) ): 1: num of beams 


% Skip over data with 0 values --> cause errors 
fem statistics 

Te eo comr (1) | == 0) ean 

else 


row_index = row_index + 1; 


% Calculate amplitude of data (dB) 
Logmcatas—— 20-7 tog l0 (abs (bmi{:) a) )); 


mean_data = mmean(log_data,N,overlap) ; 
lwad_stats_slope(row_index,:) = 
gradient (mean_data) ; 


% Detrend the data 
data = log_data - ... 
DOlyval(polyt1e (x, tog datas) -<). 


lwad_stats_std(row_index,:) = 
mstd(data,N, overlap) ; 

lwad_stats_skewness (row_index, :) 
mskewness (data,N,overlap) ; 

lwad_stats_kurtosis (row_index, :) 
mkurtosis(data,1,N,overlap) ; 


end 
end 
end 
end 


[M,N] = size(lwad_stats_slope) ; 


lwad_stats_slope((M+1),:) 
lwad_stats_slope((M+2),:) 


mean(lwad_stats_slope) ; 
std(lwad_stats_slope(1:M,:)); 


lwad_stats_std((M+1),:) = mean(lwad_stats_std) ; 
lwad_stats_std((M+2),:) = std(lwad_stats_std(1:M,:)); 


lwad_stats_skewness((M+1),:) = mean(lwad_stats_skewness) ; 


— 
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lwad_stats_skewness ((M+2),:) 
std(lwad_stats_skewness (1: 


jes 


i iW 


lwad_stats_kurtosis((M+1),:) mean (lwad_stats_kurtosis) ; 
lwad_stats_kurtosis ((M+2),:) 


std(lwad_stats_kurtosis(1:M,:)); 


Sun 


lwad_stats_slope((M+1):(M+2),:); 
lwad_stats_std((M+1):(M+2),:); 
lwad_stats_skewness((M+1):(M+2),:); 
lwad_stats_kurtosis((M+1):(M+2),:),; 


area_stats(1:2,: 
area_stats(3:4,: 
area_stats(5:6,: 
area_stats(7:8,: 
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function [y,actual_overlap] = mmean(x,N, overlap) 
SMMEAN Moving mean value. 
% MMEAN(x,N,overlap) should be used to calculate a "moving 


% average" where N is the number of points to use for each 
% interval, and OVERLAP is the ammount of overlap to uSe. 
% In some cases an asSigned overlap may not result in an 
% integer value for STEP (i.e. number of points to advance 
% after each calculation). In this case, OVERLAP will be 
% rounded to the nearest value which will result in an 
% integer value for STEP. The ACTUAL_OVERLAP can be 
% displayed if desired. 
% 
% See also MEAN, MSTD, MKURTOSIS, MSKEWNESS, LWAD. 
% Created: May 8, 1999 
% David Schalm 
% 
% Last Modified: July 26, 1999 
% David Schalm 
Te ipvalbetcfiial == 
overlap = 0; 
end 


step = round((l-overlap)*N); % Number of points to advance 


= le 

while (N+(1-1)*step) <= length(x) 
y(i) = mean(x((1+(1i-1)*step): (N+(1-1)*step))); 
5 ee ee ee 

end 


Te ieteapaqefowke —— 
actual_overlap = 1-(step/N); 
end 
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function [y,actual_overlap] = mstd(x,N, overlap) 
SMSTD Moving standard deviation 


MSTD(x,N,overlap) should be used to calculate a "moving 
std" where N is the number of points to use for each 
interval, and OVERLAP is the ammount of overlap to use. 
In some cases an assigned overlap may not result in an 
integer value for STEP (i.e. number of points to advance 
after each calculation). In this case, OVERLAP will be 
rounded to the nearest value which will result in an 
integer value for STEP. The ACTUAL_OVERLAP can be 
displayed if desired. 


dP AP ADP AP DP GO DO cP oP AO OP HP OX 


See also MEAN, MEDIAN, STD, MIN, MAX, COV, LWAD. 


% Created: May 8, 1999 
% David Schalm 
% 
% Last Modified: July 26, 1999 
% David Schalm 
joe ivetao eh —— i 2 
overlap = 0; 
end 


step = round((l-overlap)*N); % Number of points to advance 
ese 
while (N+(1i-1)*step) <= length(x) 


Vda) std(x((1+(1i-1)*step) : (N+(1-1)*step))); 
= ay cee 

end 

ise latenecpeis S92 
actual_overlap = 1-(step/N); 

end 
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function [s,actual_overlap] = mskewness(x,N, overlap) 
SMSKEWNESS Moving Skewness. 


% MSKEWNESS (x,N,overlap) should be used to calculate a 

% "moving skewness" where N is the number of points to use 
% for each interval, and OVERLAP is the ammount of overlap 
% to use. In some cases an assigned overlap may not 
result 

% in an integer value for STEP (1.e. number of points to 

% advance after each calculation). In this case, OVERLAP 
% will be rounded to the nearest value which will result 

% in an integer value for STEP. The ACTUAL_OVERLAP can be 
% displayed if desired. 

% 

% See also SKEWNESS, MMEAN, MSTD, MKURTOSIS, LWAD. 

% Created: May 8, 1999 

% David Schalm 

% 

% Last Modified: July 26, 1999 

% David Schalm 


I Pease ere ———s 
Overlap e— 0; 
end 


step = round((l-overlap)*N); % Number of points to advance 
1) = Abe 
while (N+(1-1)*step) <= length(x) 
s(i) = skewness (x((1+(1-1)*step) : (N+(1-1)*step))); 
a eee ee 
end 
Lf Nargewe. —— 


actual_overlap = 1-(step/N); 
end 
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function [k,actualloverlapl = migintosis @2,£taq,N,ovemlap) 
SMKURTOSIS Moving Kurtosis. 

% MKURTOSIS (x, flag,N,overlap) should be used to calculate 
a "Moving kurtosis" where N is the number of points to 
use for each interval, and OVERLAP is the ammount of 
overlap to use. In some cases an assigned overlap may 
not result in an integer value for STEP (i.e. number of 
points to advance after each calculation). In this 
case, OVERLAP will be rounded to the nearest value 
which will result in an integer value for STEP. The 
ACTUAL_OVERLAP can be displayed if desired. 


MKURTOSIS(x,1,N,overlap) calculates excess kurtosis. 


dP dP dP dO cH GP ADP AP DP DP DO AP 


See also KURTOSIS, MSTD, MSKEWNESS, MMEAN, LWAD. 


Created: May 8, 1999 
David Schalm 


Last Modified: July 26, 1999 
David Schalm 


dP dP dP dP dO 


tionakoin: == 2 
overlap = 0; 
end 


step = round((l-overlap)*N); % Number of points to advance 
ieee 
while (N+(1-1)*step) <= length(x) 

1p blag. —— 1 


Wii) ="kumtocis (xt (1+ (1-1)*step) - (N+ (a-1) *Step) ) ) same 
else 
k(i) = kurtosis (x((1+(1i-1) *step) : (N+(1-1)*step))); 
end 
_ = 2 + le 
end 


1£ nargout == 
actual_overlap = 1-(step/N); 
end 
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FUNGCEION eqlatout, eelonout,tgtlatout,tgtlenout)] = 
lwad_posit (data_file) 

S6LWAD POSIT 

Determine the coordinates of the DD (lat/lon) based on 

the ping time. 


oP do oP oP 


See also LWAD 


Created: May 15, 1999 
David Schalm 


Last Modified: July 26, 1999 
David Schalm 


dP oP dP dP dP 


load (data_file, ’-mat’) 

load evt31.mat 

ret (hime. —sste2num (pingqetame (1:2) )=*3600ne 
StnlmMuM ping seime (45) j-60. + Stronum (pling eime wo) - 


ale 

while evt31_time(i) <= ref_time 
1 =e 1 

end 

interp_factor = (evt3l_time(i) - ref_time) /60; 

ejlatvower= co latii) B= (cg labk(1) = -cg_lat(1-1))* 2... 
interp_ factor; 

CQloneCuUtLe segment (Co Jon (1) cq_lon(a—-iia 
interp factor; 

EqGulAabcOutLe=mtcGgt alae) = (Cotelat(1)" = tgt_latii—1))*2 
interp factor; 

EQELOMOMC = EGt elon tl) sm (tgqte lon (1) ]— ge. Ponta — 1) hoa 


interp_factor; 
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function [course,speed,ctrbrng,botdepth,deangle] = 
lwad_prm(data_file) 
$LWAD_PRM 


% This function will access the parameter file 

% ‘evt3lprm.mat’ to determine the ping parameters 
% associated with the selected data file. 

% 

% course: ships course in degrees True 

% speed: ships speed in knots 

% ctrbrng: center of ping section in degrees True 
% botdepth: bottom depth in fathoms 

% deangle: sonar depression angle in degrees 

% 

% See also LWAD 

% Created: May 11, 1999 

% David Schalm 

% 

% Last Modified: July 26, 1999 

% David Schalm 


load(data_file, ’-mat’) 

load evt3iprm.mat 

ref time = str2num(ping_time(i:2))*3600 + 
str2num(ping_time(4:5))*60 + str2num(ping_time(7:8) ); 


ele 
while timeday(1i) < ref_time 
10 ae elt 


course = course(i); 
speed = speed(1); 
CtEbEneg = Cimcomng (a) 
botdepth = botdepth (1); 
deangle = deangle(1); 
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function lwad_map (fname) 
%LWAD_MAP 
Create a map of the exercise area 
(1) A map of the LWAD 99-1 exercise will be created 
showing the following: 

(a) blue cg track history, red circle is cg 
current location 

(b) red tgt track history, blue circle is tgt 
current location 

(c) ping section is shown in yellow with dashed 
yellow center section true bearing line. 
Yellow asterisks show 10nm range. 

(d) latitude/longitude grid lines (black dashed 
lines) 

(e) Sediment map showing hal hook (black circles) 
with phi size approx 4-7, clay silt (black dots) 
with phi size approx 6-8, and sand 
(black slashed lines) with phi size approx 1. 
This map is based on data received from 
J. Fulford (NRL Stennis). 


If no filename is specified, a gui will open to allow 
the user to select a file. 


dP dP dP dP dP dP dO dO dP dP dP dO cP dP dP oP dP dP dP dP dP ce cP 


See also LWAD 


Created: May 15, 1999 
David Schalm 


Last Modified: July 26, 1999 
David Schalm 


dP oP oP dO dP 


naeeeoyeVico see — 

load data_drive 

[filename,pathname] = uigetfile([data_drive] ) ; 

fname = [pathname filename]; 

data_drive = pathname; 

save ‘e:\schalm\lwad 99-1\data_drive.mat’ data_drive 
end 


load(fname, ’-mat’ ) 


[cglat,cglon,tgtlat,tgtlon] = lwad_posit (fname) ; 

[course, speed,ctrTbrng,botdepth,deangle] = lwad_prm(fname) ; 
load evt31.mat; 

num_of_beams = size(bmf,2); 

figure (100) 
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eli 

See 007... 
‘NumberTitle’ One, =e 
‘Name’, ’LWAD 99-1 Event 31’); 


% Plot the sediment 

pic = imread(‘florida’, ’bmp’); 
sed_lon = linspace(-84.5,-83.9,600); 
sed_lat = linspace(25.4,25.9,600); 
image (sed_lon,sed_lat,pic); 

ho lid fon acdisicleennmaccrS macy 


spiGGryeO angetger Erack history 

pDiloGl—eq. Vonveq. Vat), hold on 

xlabel(’Longitude (W)’), ylabel(’Latitude (N) ’) 
title([’LWAD 99-1 Event 31’]), axis equal, axis image 
plotl-Ccqlon,calac. 7.0 } 

DICE —egraton, tcge Tat, rr.) 

DIlOE(=tagelom caoclar. bo.) 


%$ Plot a 10nm range circle, and ctr bearing line 
theta_one = ctrTbrng - num_of_beams*5/2; 
theta = theta_one:5: (theta_one+(num_of_beams-1)*5); 


theta(num_of_beams+1) = theta_one; 
for i = 1:1: (num_of_beams+1) 
Girelewlaria) = ecqliate+ V0*cosi{ theta (2) pi so0y7/ 60; 
ci rclewton( 1). — —(eqlon. — 10*sin(theta(i)*p1/180)/ +... 
(60> Cosi catat-p1/180 i) ) ; 
end 
if num_of_beams == 24 


DlLOe(eGipelemvon, circle Jat,’ y**) 
ploee@l—eqleon Circle lon({13)),... 


[cglat circle_lat(13)],’--y’) 
Garcle late—=m[ecqtlat circle= lat tik 24) colaktl; 
ecarcie lon gwl-colon circle Jom 724) @-calon|-; 
Plot (circle Nem, cipele lat,’ y ) 


else 
Plot (Circle lomacimel eat, yy") 
plot ([-cglon circle_lon((num_of_beams/2) + 1)],... 
leqlat. carclesilat((munsotebeamsc/2) etal) i. =—y a) 
Plot (circle lon, cireteslaty, 4”) 
end 
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function [cout] = lwad_contour (fname, interval, overlap) 
%SLWAD_CONTOUR 


dP dP dP dP oP oP dO cP dP dO dP dP GO GO dP oP dO dP oP CO 


JP co dO dP of 


Create a Sediment contour based on the selected ping. 


This plot is created by determining the difference 
between area averaged slope data and slope of 
individual beam data. 


If no filename, interval, or overlap are specified, a 
gui will open to allow the user to select a file. 
The interval and overlap will default to 1 sec and 
.75. Modify this file in order to change the 
defaults. 


hhesplet shows the Cg Current position (black circle), 
target current position (red circle), and exercise 


area (black dashed lines). 


Tite mivMeGtToOnmuseorecode from the Matlab ssumction 
TEOREOUr me 


See also LWAD_CONTOURF, CONTOUR, LWAD. 


Created: May 15, 1999 
David Schalm 


Last Modified: July 26, 1999 
David Schalm , 


lipemia @/ eae = = 


load data_drive 

[filename,pathname] = uigetfile([data_drive ‘*.*']); 
fname = [pathname filename]; 

data_drive = pathname; 

save ‘’e:\schalm\lwad 99-1\data_drive.mat’ data_drive 


%%%4%E% Default values %%%%% 


interval = 1; 
overlap = .75; 
spacing = interval - interval*overlap; 


sunecenvale=—1npue( What 1S thesintervali(.S 71,2): ")-; 


SOVenlap = lnouln( “Nnate 1s seheroverlan (0, .25,25;..75).: 4); 
EESEEEEEEEEEESEEEEEESEEESEEES 

else 
1f fname((length(fname)-2):length(fname)) == ‘mat’ 


filename = fname((length(fname) -23):length(fname) ) ; 


vial 


else 


Filename = fname((length(fname)-19) :length(fname) ) ; 


end 
end 


%* Load the data file 
load(fname, ’-mat’ ) 


[cglat,cglon,tgtlat,tgtlon] = lwad_posit (fname) ; 
[course, speed, ctrTbrng, botdepth,deangle] = lwad_prm(fname) ; 
if filename(5:6) == ‘SD’ | filename(5:6) == ’sd’ 
file_type = filename(5:6) ; 
else 
file_type = filename(19:20); 
end 


% Load the average statistics 

area_stats_fname = [’e:\schalm\lwad 99-1\area_stats\’ 
file_type ’_stats’ num2str(interval) ’.mat’]; 

load(area_stats_fname) 

it Mevet lap == 


area_stats = stats_0; 
elseif overlap == .25 

area_stats = stats_25; 
elseif overlap == .5 

area_stats = stats_5; 
else 

aréa_stats = stats_7/5; 


end 


% Create time\range vector 
sound_speed = .81; %nm/s 

N = round(interval*samplerate) ; 
data. start = 1; 


data_duration = length(bmf)/samplerate - (1/samplerate)... 
+ data_start; 

time = (data_start:1/samplerate:data_duration) ’; 

[mtime,ao] = mmean(time,N,overlap) ; 

mrange = sound_speed*mtime/2; 


% Create bearing vector 

num_of_beams = size(bm£,2); 

theta_one = ctrTbrng - num_of_beams*5/2 + 5/2; 
theta = theta_one:5: (theta_one+ (num_of_beams-1) *5) ; 


index Oe 
for 1 = 1:1:num_of_beams 
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ie sumneomedie: , 2) =—- 0). ond 
% skip this beam 


else 
hegedeata = 20* logl0 (abs (bmer (2,1) )- 
index = index + 1; 


thetal (index) = theta(i); 
mean_data = mmean(log_data,N,overlap) ; 
slope_data = (gradient (mean_data) ); 
diff data(index,:) = slope_data - area_stats(l1,:); 
end 
end 


cl = contourc (mrange, thetal,diff_data) ; 


ie = ele 

while i < length(cl) 
nexti = i + c1(2,i); 
Guat ieee le eerie 


% convert range/brng to lat/lon 
fOr ake) (iit Lb) 2 1 onexti 
eG newiels Kj =e— (Calon — Ssinici(27k) ora 0) =... 
Cite oe cos. (colatcsor, 1 s0jea 
emow (2) 1) s=meCOlac teecos (eli 2.) ~p1/180)* 2. 


Gre a aie) e/o 10). 
end 
2 = next, + 1; 
end 
eo = Cenc. 
figure 


Mims r= sie e4e4  —s4 25.5 252761; 


colortab = get(newplot, ‘colororder'’); 
[me newt v= "Size (Colonptab):; 


% Create the plot 
newplot; 
if ~ishold 
View (3) grid en 
Set(gqca, <lim elims@e-2) 7 ylam sims (3-4) 
end 


damit %="“sizelc,2)- 
i=. 
elle 


SOlLOro nN = "ge 


25 


while(i < limit) 
Zz level =-c(13 aie 


MPO lnteS  =—e 2 aan 

nexti = it+npoints+l; 

xdata = c(1,i+1l1:i+npoints) ; 

ydata = c(2,1i+tl:i+npoints); 

zdata = z_level + O*xdata; % Make zdata the same size 


% as xdata 


% Create the patches or lines 

cu = patch(’xXData’, {[xdata NaN], ’YData’, {ydata NaN], 
‘ZData’, [zdata NaN], ’CData’, [zdata NaN], 
‘facecolor’,’none’,’edgecolor’,’flat’,... 
‘userdata’,z level); 

= eeu); 

color _h = sieoler ne-iz level]; 

i= escka: 

end 


for, = Jelengtn (hm) 
seth je. 4oata 7.11): 
end 


view(2); 
SeriGca,” box’, On’) > 


axis equal, grid on, zoom on, hold on 
axis([-84.5 -83.9 25.4 25.9]) . 
xlabel(’Longitude (W)’), ylabel(’Latitude (N)’) 
warning off 
Ei ble@boeven= 30/¢e filename (13:16)-" > * 

Pang date "|" ping_time] ) 
caxis([-12 12]) 
colorbar 
plot (colon, collate. ko7 ) 
plot (-tqtlon, tagelat, noe 


we rlot cq box 

DIloE(|—84.27 =84.27 84.11 ,—-84.11 -84. 2/7 = 
5 61925274) 2 See Gb 25.61) eee) 

Zoom on 

brighten (1) 


if nargout > 0 


COuUne— eC: 
end 
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function [cout] = lwad_contourf (fname, interval, overlap) 
SLWAD_CONTOURF 


% 


dP dP dP dP dO dP of dP dO dP dO dP dP oP DP dO dP dO oP oP 


dP of dP dP dP 


Create a filled sediment contour based on the selected 
ano. 


This plot is created by determining the difference 
between area averaged slope data and slope of 
individual beam data. 


If no filename, interval, or overlap are specified, a 
gui will open to allow the user to select a file. 
The interval and overlap will default to 1 sec and 
.75. Modify this file in order to change the 
defaults. 


The plot shows the cg current position (black circle), 
target current position (red circle), and exercise 
area (black dashed lines). 


This function uses code from the Matlab function 
“Conteomr ft .m’ 


See also LWAD_CONTOUR, CONTOURF, LWAD. 


Created: May 15, 1999 
David Schalm 


Last Modified: July 26, 1999 
David Schalm 


ied haan — 


load data_drive 

[filename, pathname] = uigetfile([data_drive '*.*’]); 
fname = [pathname filename] ; 

data_drive = pathname; 

save ’e:\schalm\lwad 99-1\data_drive.mat’ data_drive 


$%%$%%% Default values %%%%% 


interval = 1; 
overlap = .75; 
interval = input(’What is the interval (.5,1,2): ‘'); 
SOveGlap =—inpule( Whatwas theseverlap (0,-25,.5,./5): '); 
EEEEEEEEEEEEESEESEEEESEESEES 

else 
if fname((length(fname)-2):length(fname)) == ‘mat’ 


filename = Peet Lengiehamame) 23) length ( fname) ); 


oD 


else 
Filename = fname((length(fname) -19) :length(fname) ) ; 
end 
end 


% Load the data file 

load(fname, ’-mat’ ) 

[cglat,cglon,tgtlat,tgtlon] = lwad_posit (fname) ; 
[course, speed, ctrTbrng,botdepth,deangle] = lwad_prm(fname) ; 
spacing = interval - interval*overlap; 


if filename(5:6) == 'SD’ | filename(5:6) == ‘sd’ 
file_type = filename(5:6) ; 

else 
file_type = filename(19:20) ; 

end 


% Load the average statistics 


area_stats_fname = [’e:\schalm\lwad 99-1\area_stats\’ 
file type ‘’_stats’ num2str(interval) ’.mat’]; 
load(area_stats_fname) 
ihe verlap. == 0 
area_stats = stats_0; 
elseif overlap == .25 
area_stats = stats_25; 
elseif overlap == .5 
area_stats = stats_5; 
else 
area_stats = stats_75; 
end 


%$ Create time\range vector 
sound_speed = .81; %nm/s 

N = round(interval*samplerate) ; 
data_start = 1; 


data_duration = length(bmf)/samplerate - (1/sSamplerate)... 
+ data_start; 

time = (data_start:1/samplerate:data_duration) ’; 

[mtime,ao] = mmean(time,N,overlap) ; 

mrange = sound_speed*mtime/2; 


% Create bearing vector 


num_of beams = size(bmf,2); 
theta_one = ctrTbrng - num_of_beams*5/2 + 5/2; 
theta = theta_one:5: (theta_one+(num_of_beams-1)*5) ; 
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index 0; 
Lon 1 1:1:num_of_beams 
ieee olnlonie (ee) —— 60). Sa) 
% skip beams w/ bad data 
else 


log_data = 20*logl0 (abs (bmf(:,1))); 


index = index + 1; 
thetal (index) = theta(i); 


mean_data = mmean(log_data,N,overlap) ; 
slope_data = (gradient (mean_data) ); 
diff_data(index,:) = (slope _data - area_stats(l1, 


end 
end 


diff_data = diff_data/spacing; 


ieee 

Cole 

x = mrange; 

y = thetal; 

Zo mecist tcat a 

nv = []; 

Siva 4 

ie cone) —— 1) a vay send: 
Ties Tet sees ==) 2 x= > cond: 


ima, mazl-="size(z); 
bimsre tt =o4.4 64) 2525 257,81. 


lao ee GaSe teen cae) 
Mine = mimi 2 (1) 
maxz = Max (zt ii 


ope Jbesavotwlaylania ik 
1f isempty (nv) 


GS-=coneounse ({mina maxz = ming maxziie 


else 


GS =contoure( |mainz. MaxZ “2 minz maz, 


end 


% Find the levels 
Pa. =. le 


nv = minz: % Include minz so that the contours are 


Stotally filled 
while (11 < size(CS,2)), 
Miv= [ryee > (era); 
P= i te eS (2) 11 + 1: 
end 


of 


ema} 


end 


$ Don’t fill contours below the lowest level specified in 
nv. To fill all contours, specify a value of nv lower 
% than the minimum of the surface. 
draw_min=0; 
ite mcltity (dln <= eine) 

draw_min=1; 
end 


re) 


% Get the unique levels 


Hives = eS Oise nen 2 rane( ee 
ZA, Sd (lease Tal let Glade te (tag) yt sie 
ny enya de, 


% Surround the matrix by a very low region to get closed 
% contours, and replace any NaN with low numbers as well. 


Zzz=[ repmat(NaN,1,nz+2) ; repmat(NaN,mz,1) z 
repmat(NaN,mz,1) ; repmat (NaN,1,nz+2) J; 

KK=findersnam( zz (-3.) ) 

zz (kk) =minz-l1e4* (maxz-minz) +zeros (size(kk) ); 

SO Se eee eX (ee) xX; 62 ees, ) sc nz). = 

Vivien Vee NY (Ze 7 2 eu) ya 

1£ (min(size(yy) )==l1), 

(eS Msg cCOntoOurs (552, yy , 22, av). 
else 


PeS- mse l—cCompeuircrnoeca(. 1 lsmz mz ise yu ot deena nZ |e 


pou 
end; 
if ~isempty (msg), error(msg); end 


c_ new = []; 

i=" I 

while i < length(CS) 
nexti = 1 + €S(27a6" 
ec. new( 2, 1) = CSics 2): 


% convert range/brng to lat/lon 
how Ki (2 nesta 
conew( 1) j= eollon — sin (€Si(2, 4) pi/ 1280). 
esi iy (od neos (cglat*pi/ 1809) 
cumew (2) Kime] comet +. cos (CSt2Z- ku, 120) 72]. 
Coil k 7/60; 
end 
i = Mexta. + > 
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end 
CS = c_new; 


figure 


% Find the indices of the curves in the c matrix, and get 
% the area of closed curves in order to draw patches 
% COLKGECElye 
ele le 
necurves = 0; 
r= [le 
Area=[]; 
while (11 < size(CS,2)), 
Ml Cs (2 ran). 
ncurves = ncurves + 1; 
Limeursvyes) j= 21; 
Moco lite pi) tee irst patch 
Mio] Saez eee le Tile) ye 
Area(ncurves)=sum( diff(xp).* (yp(1l:nl-1)+yp(2:nl))/2 ); 
dees ae. ee les 


end 
newplot; 
1£ -i1shoeld, 

view(2); 

Set(qca,  xlim gelams (ie jee Lim | dmencs 4 3) 
end 
% Plot patches in order of decreasing size. This makes sure 
% that all the leves get drawn, not matter if we are going 
% up a hill or down into a hole. When going down we shift 
% levels though, you can tell whether we are going up or 
% down by checking the sign of the area (Since curves are 
% oriented so that the high side is always the same side). 
% Lowest curve is largest and encloses higher data always. 


He is 
[FA, 1A]=sort (-abs (Area) ); 
if ~asstriget(gqca,;, color’ ))-, 
bo = Geb(cGea,) color > 
else 
be '= “GeE(oct, color or 
end 
if isempty (col) 
edgec = get(gcf, 'defaultsurfaceedgecolor’ ); 
else 
edgec = col; 
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end 
1f isempty (lin) 
edgestyle = get(gcf, ‘defaultpatchlinestyle’); 
else 
edgestyle = lin; 
end 
Or 7} \=eA, 
Hr—Co ioe (1 ae 
hey —CSih i (aa ee 
if (lev ~=minz draw_min ), 
PISA SEG BE (9) af jae (AE Gigi) yi 
Vp=€S. (2 aaa) ie nl); 
1f (sign(Area(jj)) ~=sign(Area(IA(1))) ), 
Kic—f and (mv—==lev7)- 
1f (kk>1+sum(nv<=minz) * (~draw_min)), 
lev=nv (kk-1) ; 
else 
lev=NaN; % missing data section 
end; : 
end; 





1f (isfinite(lev)), 
H= | Hea Cen ese, Vo, Wey, tacecolonm 2. Elet’....4 
‘edgecolor’,edgec, 'linestyle’,edgestyle) ]; 
else 
H=|H;paten(xpyyp lew tacecolor’ , bg, <<< 
‘edgecolor’,edgec, ‘linestyle’,edgestyle) J; 
end; 
end; 
end; 


numPatches = length(H); 
if numPatches>1l 
for i=l:numPatches 


set (H(i), ‘faceofisetEactor’, 0, .«.. 
‘faceoffsetbias’, (le-3)+(numPatches-i1)/... 
(numPatches-1)/30) ; 
end 
end 


% Set the plot options 
axis equal, grid on, zoom on, hold on 
aX LSi( toe 5S. o eee So |G) 
xlabel(’Longitude (W)’), ylabel(’Latitude (N)’) 
warning off 
title([’Ping:’ filename(13:16) 
‘, Mode: SD, Interval: 1, Overlap: .75’]) 
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% Uncomment the next line for a b/w plot 
$colormap (gray) 


caxis([{-12 12]) 

colorbar 

plot (-eqVlenmeciar, ko} 
plot (qtotlem, ecelar, “so ) 


% PIO cq bex 


DloOt(le34 2 eee 2 ee Cae ee 84 2 lianas 


[Zor Ole A. / PZ Oder oe Ge Ko”) 
zoom on 


LEY Marcel 0 


COUT = GEG: 
end 


ee 


SOMPOSE. 


This script file will create a composite graph of many 
separate pings. For now, it will only work for SD mode 
pEngs. 


oP of of dO oO 


de 


Location of the pings to combine 
filepath = {’e:\schalm\LWAD 99-1 Data\Tape 17\sd\'}; 


%$%%% Default values (Modify if desired) %%%% 
interval = 1; 

overlap = .75; 

spacing = interval - interval*overlap; 
ESEEEEESESEESEEEEEEEEEEEEEEEEEELEEESEEEEEEEEES 


% Load the area stats file 
area_stats_fname = [’e:\schalm\lwad 99-1\area_stats\’... 


‘sd_stats’ num2str(interval) ‘.mat’]; 


load(area_stats_fname) 


if overlap == 0 
area_stats = stats_0; 
elseif overlap == .25 
area_stats = stats_25; 
elseif overlap == .5 
area _stats = stats_5; 
elcetteoverlap ==] 475 
area_stats = stats_75; 
end 


% Get a directory of the pings to combine 
Davee (eharmi (takepart) ).: 
[fname{1:length(D)}] = deal(D.name) ; 


$ get rid of ’.’ and ’..’ filenames 
fname = fname(3:length(fname) ); 

% Select the Lat/Lon you want to look at 
% and the grid resolution 

lon, =s=s—e4. 45.201, .303 9; 

Later 25.045 204 25297 


composit_graph = zeros(length(lon),length(lat)); 
FeramMactrixe= ZeEros (Length (lon) Jengen (iat) ); 


= Pill anvehe vyalues@onL FeoMmposteaotrapm one ping at 4 time 
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Oma, =a lat) 
for m = 1:1: length(fname) 
load([char(filepath) char(fname(m))],’-mat’) 
Gisp([’Current file is ’ char(filepath) char(fname{m}) ]) 
[cglat,cglon,tgtlat,tgtion] = lwad_posit... 
({char(filepath) char(fname(m))]); 
[course, speed, ctrTbrng, botdepth,deangle] = ... 
lwad_prm([char(filepath) char(fname(m))]); 
cglion = =-cglion; 


fom lomo neex — 1. lengtn (lon) 
POG einae< = 1: Length (lat) 
%* check to see if lat/lon is in range of the cg 
dx = 60*cos(cglat*pi/180) * (lon(lon_index) -cglon) ; 
dy = 60* (lat(lat_index) - cglat); 
ranGewa=esoreiOxeZ 4 Gy 7): 


if range > 9.9 
%* do nothing because this is out of range 
% of the cg 

else 
brng = atan2 (dy, dx) ; 


% Convert bearings to true (north up) bearings 
Site ora Cie ane 
beng ==)" (beng 1307 51-90)- 


else 
bisticl = 9O-preng? Leu, pay 
end 
De) Disnee 70 
bungiee— rng i360; 
end 
brings 
beam_num = ceil (brng1/5) ; 
% Skip over data with 0 values --> cause errors 
% in statistics 
iP Pecii omer. Oeateahun) a ——— 0 ue 


do nothing 


%* Skip over beams that are +/- 45 degrees 
% of the stern 
elseif course > 075 & course < 105 & 
beam_num > 45 & beam_num < 63 
% do nothing 
elseif course > 165 & course < 195 & 
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(beam_num < 9 | beam_num > 63) 
% do nothing 
else 
N = round(interval*samplerate) ; 
log_data = 20*1l0g10 (abs (bmf(:,beam_num) ) ); 
mean_data = mmean(log_data,N,overlap) ; 


slope_data = (gradient (mean_data) ); 

diff_data = (slope _data - area_stats... 
(las) / Spacang- 

sound_speed = .81; %nm/s 


data_time = range*2/sound_speed; 

index = round((data_time-.25)/.25) ; 

if index < 1 
index = 1; 

end 

Gitf point = adits sdatatinde x). 

composit_graph(lon_index,lat_index) = 
composit_graph(lon_index, lat_index)+... 
Gatti poine: | 

ref_matrix(lon_index,lat_index) = 
ref_matrix(lon_index,lat_index) + 1; 

end 
end 
end 
end 
end 


% Divide the lat/lon matrix values by the number of 
% entries at each point 
composit_graph = composit_graph./ref_matrix; 


% Plot the composite graph 
figure 
Concours: (lon, lat composite Graph Je, hold on 


% Uncomment the next line for a b/w plot 
Scolormap (gray) 

caxis([-8 8]), 

Gated on, Co Lonbar 


xlabel (‘Longitude (W)’), ylabel(’Latitude (N)’) 
title([’Composite Plot’]) 


S Plot cq box 
plot([-84.27 -84s27 =@4.11 -84.11 -84.2%],.. 
S568 25771 25. 7s. 61 25,601 ee k-—- @ 


zOom Of 
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$LWAD_FILE 
% 
% Use this script to open and load an LWAD 99-1 file 


load data_drive 

[filename,pathname}] = uigetfile([data_drive]); 

fname = [pathname filename] ; 

data_drive = pathname; 

save ‘e:\schalm\lwad 99-1\data_drive.mat’ data_drive 


% Load the data file 
load(fname, ‘’-mat’ ) 


GS 


THIS PAGE INTENTIONALLY LEFT BLANK 
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